A LAX-WENDROFF TYPE THEOREM FOR UNSTRUCTURED 
QUASIUNIFORM GRIDS 

VOLKER ELLING 



Abstract. A well-known theorem of Lax and Wendroff states that if the se- 
quence of approximate solutions to a system of hyperbolic conservation laws 
generated by a conservative consistent numerical scheme converges boundedly 
a.e. as the mesh parameter goes to zero, then the limit is a weak solution of 
the system. Moreover, if the scheme satisfies a discrete entropy inequality as 
well, the limit is an entropy solution. The original theorem applies to uniform 
Cartesian grids; this article presents a generalization for quasiuniform grids 
(with Lipschitz-boundary cells) uniformly continuous inhomogeneous numeri- 
cal fluxes and nonlinear inhomogeneous sources. The added generality allows 
a discussion of novel applications like local time stepping, grids with moving 
vertices and conservative remapping. A counterexample demonstrates that the 
theorem is not valid for arbitrary non-quasiuniform grids. 



1. Introduction 

Consider the Cauchy problem for systems of first-order conservation laws 
d 

^=o ^y- 

u{0,x)=uq{x) (areR^), (2) 

where M++^ := (0, oo) xR'^, u : R'l+^ P (where P C M" is a bounded open 
subset of the set of physically admissible values), / = (/o, . . . , ,fd)' with /o(w) = w 

and smooth fluxes f, = (/a, ■ • • , fim)' ■ P x M'|+^ ^ K™ (i = 1, . . . , d), smooth 

source p — {pi, . . . ,p„iy : P x R'J+^ R™ and initial values uq ■.W'-^P. 

For the analysis of initial- value problems it is common to separate the time variable 
and the spatial variable(s); however, for the purposes of the Lax- Wendroff theorem 
there is no benefit in distinguishing them. For brevity of notation we collect them 
in the single vector y ~ {t,x')' £ W^^. x will be used for coordinates in R'', V 
resp. S for the {d + l)-dimensional resp. d-dimensional Hausdorff measure. 
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It is well-known that even for smooth initial values uq, there need not exist a smooth 
solution to for all yo > 0. It is necessary to extend the search to weak solutions, 
i.e. to u G Li(M'f+^;P) that satisfy 



7r^+i ^ % Jro 



0(0, dS+ I ^^^p{u{y),y)cf,{y)dV{y) 

(3) 



for all (f) S {M.'^^). Since there can be more than one weak solution, an entropy 
condition is needed to select the "physical" one: let 77 = (770, • . • , Jyd)', rji : P x 



M smooth, 770 (the entropy) strictly convex and rji, . . . ,rid {entropy fluxes) 
such that 

- — => — — (i = a= 1,...,to); (4) 

OUa ^ oup dua 

r] is called entropy /entropy flux pair. The entropy condition is 

^ d ^'^ d d ' 



which is meant to hold in the weak sense, i.e. for all nonnegative (f> G C^(R5j. ) 

/ Vo{uom{0,-)dV+ [ 4>{y)g{u{y),y)dV. 

(6) 



(Note: whether this entropy condition is sufficient to guarantee uniqueness is not 
known, except for some special cases.) 

The rest of this section is limited to the case of conservation laws without sources; 
the presence of sources, as in reactive flow, poses additional difficulties. 

The classical proof that a sequence («'') of numerical approximations to con- 
verges to an entropy solution proceeds as follows: the properties of the numerical 
scheme (e.g. a monotone conservative scheme with consistent homogeneous (i.e. y- 
independent) fluxes on a uniform Cartesian grid, see IHHL76| . [HMSOl and ET8OI) 
guarantee that (u'*) is bounded in L°° and TV (the space of functions with bounded 
variation in the sense of Tonelli-Cesari) . This implies that some subsequences con- 
verge pointwise almost everywhere (in the presence of uniform L°° boundedness 
equivalent to Lj^^ convergence); the Lax-Wendroff theorem (see jLW60| ) proves 
that the limits of these subsequences are indeed weak solutions of Moreover, 
if the satisfy discrete entropy inequalities, then the limits must be entropy solu- 
tions of (Q). Whenever entropy solutions are unique, the entire sequence {u'^) must 
converge to the entropy solution. Positive uniqueness results are available in special 
cases (see |Kru70| for scalar conservation laws in multiple dimensions or |IjL97j for 
ID system entropy solutions with small total variation and some other restrictions), 
but see |E1103| for a possible counterexample for 2D Euler system solutions. 
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In many cases, Lj^^ precompactness is difBcult to prove — and might be false — , 
e.g. for unstructured grids or higher-order schemes for scalar conservation laws, not 
to mention schemes for systems of conservation laws. (For this reason, techniques 
based on measure-valued solutions which require L°° boundedness, but not Lj^^. 
precompactness, have been developed and successfully applied to the scalar case 
in j(]L91j . |('L93j : see also |Noe95j for irregular grids.) However, |(X]L94j have 
generalized earlier work by |Kuz75j (see also |San83j ) to a large class of unstructured 
grids. They prove convergence of order i to the entropy solution, for monotone 
numerical fluxes with antidiffusive modifications; the modifications allow for higher 
order in regions where the entropy solution is smooth. Although these Kuznetsov- 
type proofs yield convergence without resort to the Lax-Wendroff theorem, they 
demonstrate that the preconditions of the Lax-Wendroff theorem are satisfied more 
often than previously thought. 

More importantly, while rigorous proofs of convergence are limited to special cases, 
observing actual output of good numerical schemes suggests that bounded a.e. 
convergence is rather common, even for important systems like compressible gas 
dynamics. In this sense, the Lax-Wendroff theorem has important heuristic value: 
it guarantees that the limit, if there is one, is a weak solution; moreover, in the 
presence of a discrete entropy condition, it guarantees that the limit is an entropy 
solution. Finally, the Lax-Wendroff theorem serves as a theoretical motivation 
for focusing on conservative schemes with consistent fluxes (however, occasionally 
nonconservative schemes are used in practice). 

While the Kuznetsov-type proof in |CCL94| applies to a large class of unstructured 
meshes, it relies strongly on special properties of scalar conservation laws; the same 
holds for techniques based on weak convergence and measure-valued solutions. It 
seems that only the Lax-Wendroff theorem provides at least a partial result for 
systems. 

The original Lax-Wendroff theorem requires a ID uniform Cartesian grid, contin- 
uous fluxes, LJq^ precompactness and boundedness. LeVeque | LeV92| Section 
12.4 simplifies the proof, at the cost of requiring locally Lipschitz-continuous nu- 
merical fluxes and TV boundedness. |KRW96j present a proof for 2D polygonal 
meshes, locally Lipschitz-continuous numerical fluxes, L°° boundedness, Lj^^ pre- 
compactness and an explicit CFL condition; however, their assumptions (2.3) and 
(2.4) about the mesh seem restrictive. More general triangular meshes arc covered 
by Proposition 4.4.1 in |GR96| : it might be possible to extend their proof technique 
to polygonal meshes. With straightforward modifications to statement and proof, 
all of these results and proofs apply to an arbitrary number of dimensions; however, 
none of them seem to generalize into other directions easily. This article considers 
a quasiuniform mesh with no other geometric restrictions, (uniformly) continuous 
inhomogeneous numerical fiuxes and nonlinear inhomogeneous source terms. 

Only the Cauchy problem is discussed; boundary conditions pose many open re- 
search problems, both theoretically and numerically. Even in benign cases where 
the flux is completely prescribed and independent of the solution near the bound- 
ary (as in supersonic inflow), one needs to make additional assumptions about the 
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convergence of the numerical solution near the boundary which are not implied by 
mere boundedly a.e. convergence. 

Scction[21introduccs the grids, numerical fluxes and numerical sources and the con- 
ditions imposed on them; this rather abstract framework is illustrated by a simple 
example in Section 12.61 Section |31 contains statement and proof of the generalized 
Lax-Wendroff theorem (Theorem ^ . Section 21 provides a counterexample that 
explains why Theorem ^ does not always hold for non-quasiuniform grids. The 
newfound generality enables theoretical discussion of some numerical techniques 
and applications in Sectional 



2. Notation and assumptions 



2.1. Landau symbols. Two sequences of grids will be used: unstructured grids 
with parameter h, and uniform Cartesian grids with parameter H. An expression A 
is said to be 0{B) {B some other expression) if there is some constant c, independent 
of 

C, N, F, h, H, e, p, k, w, 

so that 

A<cB 

as long as h, H € (0, 1] and as long as 

^^=i/<2 

A is said to be n{B) if B is 0{A). 

We say that an expression is Op^e(l) if, for any fixed values of p, e > (and h := pH) 
it converges to as iJ | 0. 



2.2. Grids. For any ft, > 0, let C'^ be a system of closed subsets (called cells) of 
Wi^^ with pairwise disjoint interiors so that 

U c = ]rp". 

We require the cells to have Lipschitz boundaries; this is more than weak enough 
for aU conceivable numerical meshes. For C,N e C'^ with S{C HN) > 0, let C^N 
denote the ordered pair {C,N); depending on the context it will refer to C n iV 
instead. The unit normal nc^niy) in each point y & C C\ N is fixed as pointing 
into N . The C—^N are called interior faces] the other class of faces consists of 
initial faces C— >9, d^C (where C fl ({0} x W^) ^ 0). d^C will sometimes refer to 
C n ({0} X W^), with unit normal (1,0,..., 0)'; C^d wiU refer to the same surface 
with opposite unit normal. Define C'^ := C'' U {&}. Let be the set of interior 
faces, J-^ the set of all faces. 

To "define" the mesh parameter ft, require 

diamC<ft (CeC'*); (8) 
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this implies V{C) < h"^^^. On the other hand, the mesh must be quasiuniform in 
the following sense: 

y(c) = n{h'^+'^) (c € c''). (9) 

Moreover, the cell surface measure must be controlled: 

S{dC) ^ 0{h'^) (C e C^). (10) 

Let B^' be the cr-algebra generated by C'* over R'J+^ the elements of M{B^'\ P) (the 
space of ;B''-Borel-mcasurable maps on into P) are called grid functions. For 

w'' e M (i3''; P), let e P denote the constant value of m'' on int C (C G C^). 



2.3. Numerical fluxes. Over every face F G there is a numerical flux Ep : 
M{B^) M. (Note; the following definitions make sense for numerical entropy 
fluxes. The usual numerical fluxes can be reduced to this case; see Section [2.61 ^ 
The following requirements are imposed on numerical fluxes over interior faces. 

(1) Consistency: For w P, let w be the constant grid function with value w 
(i.e. wc = w for all C G C''). We require 

Ec-,n{w)=1 r}{w,y) ■ nc^N{y)dS{y). (11) 
Jc^N 

(2) Uniform continuity: There is a function 5e : (0, cxd) — > (0,cxd) so that 

V e > 0, /i > 0, F e j^'', wG K", w'' G i°°(B''; P) : 

u)||ioo(H..p) <5B(e) |i;j^(u;'')-£;;^(u))| <e/i'^ (12) 

(again, w denotes the constant grid function with value w). 

(3) Uniform houndedness: for any^ sequence (u;'')/i>o, 

Ef{w'') = 0{h'^) {F G .F''). (13) 

(4) Conservativeness: for all w'' G M{B^;P), C^N G J^'', 

Ec^n{w'') = -En^c{w''). (14) 

(5) Bounded stencil: define the stencil oi F E as 

stni^ := {C G & : w'' ^ Ef{w^) not constant in Wq} 
and require 

sup d{C,F)^0{h). (15) 

Cestn_F 

For initial faces, we impose the numerical initial condition 
Ed^c{w'') ^ -Ec^aiw'') ^ f Muoix),0,x)dSix) V«;'' G A/(B^ R™). 

Jd^C 



(16) 



^To apply this to common cases, restrict P to be bounded (scalar case), bounded away from 
vacuum (gas dynamics) etc. 
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2.4. Numerical sources. The source terms in © are approximated by numerical 
sources: fmictions Gc '■ M{B^; P) — > R for each ceU C e C'\ The numerical sources 
must satisfy the following conditions (that are very similar to the ones for numerical 
fluxes) : 

(1) Consistency: 

yweP, CeC'': Gc{w) = / g{w, y)dV{y) (17) 

Jc 

(where w is the constant grid function with value w £ M™). 

(2) Uniform continuity: there is a function So : (0, oo) — > (0, oo) so that 

V e > 0, h>0, Ce C'\ weP, w'' G L°^{B''; P) : 

Ik" - M\l^(b'^;P) < Scie) \Gc{w'') - Gc{w)\ < eh^+K (18) 

(3) Uniform houndedness: for any sequence (u''')/,,>o, 

Gciw'') = 0{h'^+^) {CeC''). (19) 

(4) Bounded stencil: define the stencil of C G C'' as 

stnC := {C G C'' : w'' ^ Gc{w^) not constant in w^,} 
and require 

sup d{G,C')^0{h). (20) 

C'GstnC 

2.5. Result. 

Theorem 1. // a sequence {u^)h>o of grid functions satisfies the discrete scalar 
inequalities 

J2 Ec^n{u^) < Gc{u^) (V/i > 0, C G C^) (21) 
and converges almost everywhere to u, then u satisfies 

2.6. An example. For illustration, the Lax-Friedrichs scheme (for a system ut + 
f{u)x = with initial condition u{0,x) = uo{x) on a uniform ID grid with cell 
size h and uniform time steps Xh {0 < X < 1 constant)) is fit into the abstract 
framework in the previous sections. For h > 0, let C'* contain the cells C", 

Cj := [jh, {j + l)h] X [nXh, {n + l)Xh] {j e Z, n e No). (22) 

Numerical fluxes: for j G Z, ?i G Nq, 

Fg^coiw'') := / uo{jh + y) dy, (23) 







/ f{w'l^,n) + f{wK ) wK -w^^X 
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Numerical sources: all = 0. It is easy to check the numerical fluxes satisfy all 
conditions, in particular consistency. The numerical solutions m'' G Af(i3'^;]R™) are 
defined by 

p{j+i)h 

Cl := / uo{x)dS{x), (26) 

FN^ciu'') = {yh> 0, VC e C'') (27) 

which is exactly the literature definition, using different notation. 

(|27|l is a system of R™-valued equations for but it can obviously be converted 
into 2 • m systems of scalar inequalities of the type (I21() : Theorem ^ applied to each 
of them separately implies |(2J, i.e. that u is a weak solution. 

In a similar fashion, it can be verified that the limit is an entropy solution. For 
Burgers equation (/(u) = ^u^), it is sufficient to prove the entropy inequality for 
the Kruzkov family of entropies and entropy fluxes, 

Va{u) = \u- a|, ?7i(u) = sgn(u - a)f{u), (28) 

where a £ M is the family parameter. The numerical entropy fluxes 

Eq_co{u'') -.^ r^oiuoijh + y)) dy, (29) 
Jo 

i^c-^c-^^l"') -/*^?o(^&.0, (30) 

J J J 

Ec^^c^^^ {u'^) XhiFc^^c^^^ {u'^ V c) - %^c;^^ ("" A c)) (31) 

(see |UM80| ') are consistent with 77 and satisfy the discrete entropy inequality (|21|l . 
so Theorem n asserts that u is an entropy solution (in the sense ©). 



3. Proof of Theorem [I] 

The proof is based on two essential ideas: first, the original proof in )LW60j uses 
summation by parts (in analogy to the integration by parts used to derive the 
concept of weak solution); this requires a Cartesian grid. This obstacle is bypassed 
by approximating cubes with sidelength i7 in a uniform Cartesian grid by cells in 
an unstructured grid with parameter h; see Figure^ Summation by parts is carried 
out for these cubes. 

Second, since (w'') converges in Ll^^{SJ^^), u and w'' will be "close" and "nearly 
constant" in a suitable neighbourhood of "almost all" cubes (for ft, J, 0), so the con- 
tinuity and consistency properties of numerical fluxes and sources can be exploited. 
For the "few" remaining "bad" cubes, one can use uniform boundedness of numer- 
ical fluxes and sources. The proof will be completed by first fixing a sufficiently 
small ratio p to minimize geometric errors and then choosing a sufficiently small 
H > to control integral errors. 
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3.1. Cubes. Let e'*' G U^^^ (with i e {0,...,d}) be the standard basis vectors, 
with ith component = 1, aU other components = 0. Omitting the parameter H for 
readabihty, define the closed cubes 

d 

Ik:=H- Y[[h, h + 1] (fc e No X Z'^) 



i=0 



with faces 



dii^ -^hA + 1] X + 1} X n [fcj' ^-^ + 1] 

\i=0 3=i+l 

(i-l d 

note that the interiors of the Ik are pairwise disjoint and that 



U I^-^ (32) 



moreover 



a/fc = U (a/r U 9/^.+) (fc e No X Z^). (33) 



i=0 



3.2. Cube approximation. For given H,h > {h < H/2) and fc G No x Z'', we 
select an approximation Ik C C'' to Ik by requiring that 

C n 4 7^ VC e 4 (34) 

and that the sets Ik form a partition of . (These conditions need not determine 
Ik uniquely; the particular choice is not important. (|34|l admits the existence of 
such Ik because the Ik cover M^J,"*'^.) 

We also define (for i e {0, . . . , d}, fc G No x Z'*) 

dik := {C^d eT'':Ce Ik} U {C^N e T'' : C e Ik, N ^ Ik}. (35) 

gp± UC^N eT'^ -.C elk, N e 4±,,.,}, fc ± e« G No x Z-^ 
'= ■ 1{C^9g.F" :CG/(o^fe,^...,fc,)}, else 



Note that 



however 



U(5/ru5/^+)==94 

is not true in general (see Lemma 13 because some faces belong to "corners" rather 
than sides of the approximated cubes (see Figure ^1. 
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3.3. Geometric estimates. For r > and A C M'Jf^^, define the closed neigh- 
bourhoods 



Lemma 1. 



[jCcQ^ih), hcQ^iljC), (37) 
Vi\J C)= H'^+^ + 0{pH'^+^). (38) 



Moreover 



U FcQ^idlf). (39) 



Proof. These are immediate consequences of ^ and of (|34() . (|35ll resp. (|36|l . □ 
Lemma 2. A C M'|+^ meeis < ^'^fj^'i^^^ ceZfe in C'\ 

Proof. diamC < h, so C f) A ^ (d imphes C C Q^(A). However, ^ imphes that 
Qhi^) cannot contain more than V {Q i^{A))h^'''^^^'> cehs. □ 



Corollary 1. 



sup #stni^ 0(1); (40) 

sup #{7VeC'' :CniV^0} = O(l); (41) 

ceC^ 



for all fc e No X Z"* and i e {0, . . . ,d}, 



sup = ©(p-''-^), (42) 

sup #dil' = 0(p"''). (43) 

fcgNoxZ'',ie{o,...,d},se{+,-} 

Proof. (|40|l : by (|15|l . there is a constant c (independent of h) so that 

CestnF 

for all faces F. Since diami^ < /i, Qch{P) is contained in some closed ball with 
diameter 0{h), hence volume 0(/i'*+^). At most 0(1) cells fit into this ball, so 
#stnF = 0(1). 

(|41|l : this is immediate from Lemma |21 
(E^ : l|T7|) implies 

Q^( y C)C02;,(/fe), 
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Figure 1 . A cube Ik (checkerboard-shaded grid) is approximated 
by a cluster (thick boundary) of grid cells (thin triangles). 



hence {p < ^) 

V(Q,{[j C)) = 0{H''+'). 

Hence by Lemnia|21 at most 0( ^^+1 ) = 0(p~''~^) cells meet Uce/fc ^■ 
(gSJ: from derive 

Qhi U F) ^Q2hidli^), 

so 

Vm U F))^0{hH'')^0{pH'+'). 

Hence Lemma|2shows that at most Q( ^^+i ) — 0{p^'^) cells meet {Jp^gjH F; by 
(|1T|| each has 0(1) faces. □ 

The following lemma states that, for small p, "most" of dik is composed of the 
dil^ (i = 0, . . . , d), i.e. we can ignore the "corners" of Ik- 

Lemma 3. Define the "half- cylinders" 
Zf := {y e : y, ^ H{k, + i), yj e H ■ {k, + p, fc, + 1 - p) {j ^ z)} (44) 



A LAX-WENDROFF TYPE THEOREM 

















1 h 1 






i 








H 
























Ik 






















R 
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Figure 2. V{R) = 0{h'^H'^-^) = 0{p^H'^+^) 



(see Figure\^. Then 



Moreover, 



^ S{F)^0{p)H^ 



Proof. (See Figure EJ) Let F = C-^N G dllf for some i £ {0, . . . , d}, s G {+, 
(or F = C^d). By whenever F ^ Z^^ ^ F, then F (and hence C) meets 

d 

i?:=Q,(a/fe)-|J(Z^+UZr). 

i=0 

However, it is easy to verify that 

V(Q^{R)) = Oih'H''-') = 0(p'H^+'). 
By Lemma 13 at most 



oip'-") 



cells can meet R. By (|10|l . their total surface measure is 

= 0{p^-'^)0{h'^) = 0{pH'^); 

this implies H45|l. 
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Regarding H46(l : whenever C—fN G dik (with C £ /fc) is not contained in any <9/^ , 
then C—^N belongs to a "corner", i.e. G Ik+m for some m S Z''"*'^ with |m|oo = 1, 
\m\i > 2 (because of p < ^ and diam C, diam < h). This means C^N C i?; 
states that these faces may be ignored at a cost of surface measure of 0{pH'^). 

C^d e dil~ for some k e {0} x Z'', so initial faces do not contribute to the sum 
in (gBJ. □ 



The numerical fluxes over lJ_Fga7'.='= ^ ^'"^^ pieced together to approximate the 
exact flux over dll^. This requires the following geometric estimate. 

Liemma 4. For all fc £ No x Z'', i £ {0, . . . , d}, 



n dS 



0{p)H'. 



(48) 



Proof. (Sec Figure 13) 



n dS 



The first summand in H49|) equals 



ndS + 0{p)H'^ 



(49) 







/ ) 



up to 



ndS — I ndS 
Hi Jr2 



where Ri , i?2 are contained in 



dS 



(50) 



dzirnQ.idir), 

a set with total surface measure 0{p)H'^] hence the two integrals in H50|l are 0{p)H'^ 
too. □ 



3.4. Completion of the proof. Consider an arbitrary nonnegative test function 
4> G C^{W^^). Since has compact support, it is sufficient to consider the finite 
subsets 

K := {/c e No X Z^ : 3£ e Z'^+i : |^| < 1, supp n h+e + 0} 
of cube indices. Note that 

#A' = 0(p-('^+i)). 
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Qhidin 



Dpedir F - Z] 



Fed Ik 




Figure 3. The surfaces Ri, i?2, Uj^ea/r ^"^r^ and dll7 -Z]: 
have measure 0{p)H'^ and can be neglected. 



Define 



stnfc:= IJ stnFU |J stnC, 
stn A' := (J stnfc. 



keK 



Lemma 5. For all w G Lj„^(R^J+^), 



E 



Uceatnfc 



\wiy)\dV{y)^0{l) 



keK Uceatn k 



UcPntn K ^ 



(51) 
(52) 

\w{y)\dV{y) (53) 



Proof. Due to bounded stencils (|15|l resp. H20() . bounded diameters (jSJ and p < \ 
(see 0), 



sup e : C e stnfc} = 0(1). 



Hence 



E / = E *{keK:Cestnk} f \w{y)\dV{y) 



CSstn/f 



<0(1) ^ / \w{y)\dV{y). 



□ 



14 



V. ELLING 



We need to show that u, are "almost constant" and "close" on "most" cubes. 
We introduce a new parameter e > 0. Again omitting e, p, H > from the symbols 
for readability, define 

(54) 

Ik 

^ -fe|>min{^B(e),<5G(e)}}, (55) 







Bi 


= {keK 


B2 


= {keK 


B 


= Bi U ^2 


G 


^ K ~ B; 



V{y e h ■ \u{y) -Uk\> min{5£;(e), (5G(e)}} 



> e}, (56) 



Vih) 

(57) 



B contains the "bad" , G the "good" cube indices. 
Lemma 6. For any choice ofe,p>0 



T) = Op,.(l)- (58) 



Proof. Bi and B2 are treated separately. First B2. let w G C°°(Ui^^) be arbitrary 
(it will approximate w); define 

Wk = V{hr^ [ w{y)dV{y). 
•Ilk 

Then 



I \u{y)^Uk\dV{y) 
< ^ / \u{y)-Uk\dV{y) 

he 

<E/ Hy)~My)\dV{y)+Y, f Hy) ^Wk\dV{y) + 



= Oi\\u - u;|Ui(u^^^ + \\Dwh^^^^^^ j,)H) (59) 
because 

YV{h)\uk-Wk\ = Y j u{y)~w{y)dV{y) 

k£K keK •'li' 

-Y l"(y)-w(y)|dT/(y) < ||M-uj||ii(Li^^^^,j,) 

kGK •'Ik 



and because 



\wk - w{y)\ = 0{H\\Dw\\L^f<^^^^i^)) 
for y ^ Ik, by smoothness of w. For each k G B2, 

\u{y)~Uk\dViy) > emin{SEie),SGie)}ViIk), (60) 
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by definition (|56|l of B2, so H59(l implies tliat 

,^ , llii — wll 7- 1/1 I J, \ + \\Dw\\t ^(\ \ T,\H 

#B, ■ iJ'^+i = V V{h) ^ " . , (U.e.-^^) )^ (61) 

emin{(5B(e),(5G(e)} 

By first clioosing w G C°°(]R!^'*'"'^) witli sufficiently small ||u — w||li(Ua. K^k) ^^"^ 
then choosing an upper bound for H , the right-hand side can be made arbitrarily 
small. 

Regarding Bi, 



E / \u\y)~u^\dV{y) 



< E / \u\y)-uu\dV{y) 

<E/ \u\y)-uiy)\dViy)+J2 [ \u{y) - w{y)\dV{y) 

E/ \w{y)~Wk\dV{y)+y^ \wk-Uk\dV{y) 



1* 0(||w''-M||ii(u^^^,„^c) + \\u ~ wU^ii^^^^^^^c) + H\\Dw\\L^.(x}ce.^r.KC])- 

(62) 

For each k £ Bi, 

W\y)--Uk\dV{y) = f}(/+ii/'^+imin{(5B(e),<^G(e)}) = V(/fe) min{(5B(e), ,5G(e)}), 



Ucestiifc <^ 

so (inill (with ®) shows 



#Bi • H'^+^ = O 



p''+imin{5£(e),(5G(e)} 



By first choosing a suitable w and then choosing an upper bound for H , the right- 
hand side can be made arbitrarily small. □ 

Now we can finish the proof of Theorem ^ Sum 121|l over C £ Ik, multiply with 
(l){Hk) (> 0) and sum over k G K: 

E <l^{Hk) E Gc(«') > E <^(^^-) E ^^("") (63) 

(here we used the conservation property H14|l to eliminate F ^ Ufce-ff ^ ~ 

C^N with C, iV e 7fe for the same k). 

Collecting the terms on the right-hand side of 1)63(1 where F is an initial face yields 

E E Ec^o{u'^)HHk) 
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which equals 

- / rjo{uoix))(biO,x)dS{x)+0{H), (64) 
using H16() and smoothness plus compact support of 4>. 

There are at most two terms per interior face in they can be written (using 
(O) as 

Ec^N{u'')i(j){Hkc) - ^{HkN)) (65) 
where C elkc, N e 4^ (l^c - fcwioo = 1)- 

Note that the cj) difference is 0{H) (by smoothness of 0), and that for each k ^ K, 
(|46|l allows to drop all terms for interior faces that do not belong to some 9/^^, 
at the cost of terms of size 0{p)H^ ■ 0{H) = 0{p)H'^^^ per k, hence 0{p) overall. 
Here, it is important that Ef{u^) = 0{h'^) (by 

Moreover, by Lemma |S1 the number of bad cubes is Op 

^(l)i/-(rf+i)^ and due to 

uniform boundedness ifT^ . ifTUI) and 

^ \Ec^j,{u'^)\=0{H^). 

Since the (j> difference supplies an extra 0{H), the sum of terms in H65() where kc 
or kN is in B is Op,e(l)i?"(''+^) • 0{H'^) ■ 0{H) = Op,e(l). Hence the interior face 
part of is 

d 

= EE E Ep{u %<l,{H{k + eW)) - 0(gfc))^ + Q(p) + Op,.(l) 

i=0 k&G pfzQp- 

= EE E Ep{u^')H-' I ^{y)dV{y) + 0{p)+OpM) 



(using smoothness and compact support of cj)) 

= EE E EFiTl,)H-^' f ^(y)dViy) + 0{p + e)+OpA^) 

i=OkeGF<^dp- ''''' ^' 



(using the definition H57() of G in uniform continuity H12|l . combined with (|43() ') 

= J2J2 [ '?(^'^' •) ■ " d^^^' I ir^vWiy) + 0(p + 6) + op,,(i) 
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(using consistency CU) 

= E E v{^k,Hk) ■ [ ndSH-' [ ^{y)dV{y) + 0{p + e) + Op,,(l) 



i=0 k£G •^^FeoP 



(using smoothness in y of rj (note (|^ .l(7|l) with and itTUI) ') 



ED 



^^r;,(^I,,m-) / |idT/(2;) + 0(e + p) + 0p,,(l) 

EE / |^(y)^.Ky),yWy) + o(e + p) + op,e(i) 
.-n "y* 



1=0 keG 



using the definition ifsTjl of G and smoothness in y and u of r/) 
"^-EE / |^(2/)'7.Ky),yW2/) + 0(6 + p) + o,,,(l) 

= - I E l^(y)'7K^(2/)' y)'^^(y) + + p) + «p.^(i)- 

JR^+i i=0 y^ 



(66) 



It remains to treat the left-hand side of 163() . As for the flux integrals, we may 
omit the terms for fc G at a cost of Op,, (l)i7~('*+i) • 0{H<^+^) = Op,e(l) (from 
rcsp. (O), hence: 

^0(iffc) ^ Gc(«'') 
= ^0(i7fc)^ Gc(u'') + Op,,(l) 

^ (/.(i/fc) ^ Gc(^fc) + 0(e) + Op,,(l) 



ED 



V 0(iffc) / g(?Ifc, y)dF(y) + 0(e) + 0^,^(1) 
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^Y.cj){Hk)vi U c]giuk.Hk) + Oie)+OpAl) 



^ / <j>{y)dV{v)g{uu,Hk) + 0{e + p) + Op,,{l) 



= Y. I <l>{y)9{n{y),y)dV{y)+0{e + p) + 0p.,{l) 
•^E / 4>{y)9Hy),y)dViy) + 0{e + p) + OpM) 



mgHy). y)dV{y) + 0{e + p) + Op,,(l) (67) 

d+1 



3.5. Conclusion. Combining (|64|1 . (|66|) and (I67fl . we get 

0(y)5(^^(2/), y)rf^(y) + o(e + p) + o,,,(i) 



+ 



> - 



J2^.^y)^^i<y))dyiy)- / %K(x))0(o,x)d5(x) 



d+i — ' dv' 

i=0 " 



Now, we can first make the 0(e + p) term arbitrarily small by choosing appropriate 
e and p; after that, the o term can be made arbitrarily small as well by picking H. 
Therefore, u satisfies the proof is complete. 

4. A COUNTEREXAMPLE FOR NON-QUASIUNIFORM GRIDS 

While most assumptions in this paper are rather weak, an important exception is 
quasiuniformity (in the sense of lO)- Unfortunately, there is a strong counterex- 
ample to Theorem n for non-quasiuniform grids (see Figure 01): 

Example 1 (Staggered Lax-Friedrichs). Consider the trivial problem 

ut = 0, uo = X[o,i]- 
We discretize it naively with the staggered Lax-Friedrichs scheme, for flux / = 0, 
taking At = where /i > is the spatial cell size: let contain the cells 

:= [2nh, {2n + l)h] x [jh, {j + 
O^" := [(271 + l)h, {2n + 2)h] x [{j + ^-Jh, {j + ^)h] {n € No, j € 1). 
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i, n 



Figure 4. Staggered Lax-Friedrichs with At = , Ax = h: ex- 
cessive refinement in t direction causes oversmoothing and conver- 
gence to a non-solution. 

Define 

Fo^^E-^^i^'') - Fo--.E"+^ = ^"o" (" e No, J e Z), 

l'{2j+l)h 

Fq^eo{u'')= uoix)dx O'eZ), 

' J2jh 



Fej^e^^, {u'^) = 0, %^0^+, (w^) = (n e No, j £ Z), 

u^o := h^^ / (j + l)/i7io(a;)(iS'(x). 
Jjh 

It is easy to check that the numerical fluxes are consistent and satisfy the initial 
condition; all other requirements are satisfied as well. However, the (uniquely 
determined) solutions {u'^)h>o at a fixed time t approximate G{^,-) * X[o,i] (as 
/i I 0), where G is the heat kernel (this is easy to prove by considering two steps of 
the numerical scheme: 



Mr 



4 4/i 

this is a well-known finite difference scheme for ut = j^u^x, so standard theory 
applies). However, ^ — > cx), so the converge in L^^^ to 0, not to the actual 
solution X[o,i]- 



The example is so "economical" that it rules out any conceivable relaxation of the 
quasiuniformity requirement: the cells are identical rectangles (in particular they 
are convex and have well-behaved surfaces) with sides parallel to the coordinate 
axes, numerical fluxes through a face depend only on one of the directly adjacent 
cells, the numerical solutions are uniquely defined. The only "violation" is that 
cells become small in time direction faster than in space direction. 
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The example reflects a problem that does not appear in discretizations of semidis- 
crete schemes with artificial viscosity; while the viscosity coefficient of Lax-Friedrichs 
type schemes is a/i^/At (for some constant a), it is typically ah for semidiscrete 
methods, so no harm can be done by choosing At too small. 

Note that |Noe95j presents a related counterexample, namely the Lax-Friedrichs 
scheme on a uniform Cartesian (non-staggered) grid with At/h j Oas/i j 0. 
In contrast to staggered Lax-Friedrichs, this counterexample does not serve our 
purposes: for such a scheme, the flux between two intervals in a time step is not 
0{At), as required by uniform boundedness (fT!^ . but 0(1). 

The restriction to quasiuniform grids is a serious one; results for non-quasiuniform 
grids arc highly desirable because such grid sequences are produced by adaptive re- 
finement and/or adaptive time integration. However, Theorem Q probably remains 
true in a special case: tensor products of quasiuniform grids, i.e. 

s s 
a—1 a—1 

where each is a grid of the type defined in Section 12.21 (for a = 1, the grid has 
to cover R^^ , for a > 1 the grid covers R''° ) . The case of semidiscrete schemes can 
be reduced to the tensor grid case. These questions will be explored in forthcoming 
work (see |Ellaj ) . 

For arbitrary non-quasiuniform grids, on the other hand, it is necessary to impose 
stronger conditions on the numerical scheme. For example, one could study the 
error estimators that are used by adaptive schemes to determine where to refine 
the grid or decrease the time step, in order to derive additional smoothness or con- 
vergence information about {u^)h>o- If sufficiently weak assumptions can be made 
about a large class of error estimators, it might be possible to derive a Lax-Wendroff 
type result for non-quasiuniform grids that is general enough to be interesting. 



5. Novel applications 

5.1. Local time stepping. In order to resolve shocks or contact discontinuities 
well, it is necessary to refine the grid near them. The time step is limited by the CFL 
condition in small cells near these discontinuities and might be unnecessarily small 
for other parts of the domain. For this reason, it can be efficient to use different 
time steps in different regions; some schemes in this spirit have been proposed in 
|OS83j or |E1100j Chapter 4. Theorem is not limited to spatially unstructured 
grids; grids can be unstructured in space-time, as long as they are quasiuniform in 
the sense of ®. 

5.2. Moving vertices. Another use for the generalized Lax-Wendroff theorem is 
the large class of numerical methods with unsteady grids (see |E1100| Section 2.1 
for adaptation of classical approximate Riemann solvers to this case). These are 
important because some applications have moving domain boundaries, e.g. due to 
wing flutter or rotating turbine blades. Moreover, it is often natural or (for high 
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Mach number supersonic flow) more efficient to use Lagrangian methods (grid ver- 
tices move along with the fluid). To accomodate these methods is straightforward: 
instead of a tensor product of time axis partition and fixed spatial grid, a grid with 
moving cells is used; faces are no longer cither perpendicular or parallel to the time 
axes. (However, whether Theorem is applicable depends on other details of the 
scheme as well.) 

5.3. Conservative remapping. In numerical computations, it is sometimes nec- 
essary to change grids (remapping), for example because the old grid has developed 
singularities (especially common for Lagrangian schemes when there is strong vor- 
ticity in the flow field). See |Duk84l IHKWn Kiram VU^^ IDHOOj for remapping 
algorithms and applications. The remapping step should be conservative, for the 
same reasons that numerical schemes are conservative, and conservative quantities 
should not be "transported" during the remapping step more than necessary. A 
simple way to achieve this is to set 

UN = V{N)~^ V{0 n N)uo 
o 

where is a cell in the new grid, O runs over the cells in the old grid and uq, mat 
are densities of conserved quantities in each. (However, to achieve higher orders 
of accuracy, it might be necessary to compute polynomial or spline reconstructions 
v{x) from the cell averages uq and to set 

UN '■= J v{x)dV{x); 
Jn 

the previous scheme corresponds to the obvious piecewise constant reconstruction.) 

It is not clear whether remapping can prevent an otherwise fine numerical scheme 
from converging to the entropy solution. However, Theorem^ can be applied to 
answer this question for conservative remapping: the remapping step is interpreted 
as an extra hyperplane of faces (perpendicular to the time axis), with numerical 
fluxes deflned depending on the remapping algorithm. The following requirements 
are weak enough to cover most existing methods: 

(1) Consistency: whenever all uq' are constant = w, the reconstruction v 
should be constant = w in each cell O. 

(2) Continuity: the reconstruction map u'^ ^ v should be continuous on the 
"diagonal" of constant grid functions in the L°° L°° topology. 

(3) Boundcdness: the reconstruction map should be uniformly bounded; more 
precisely, for any M there should be a constant c so that 

sup|uo|<A/ [ \v{x)\dV{x) <cV{0). 

o Jo 

(4) Locality: the values of v over a cell O should depend only on cell averages 
uo' for d{0' , O) < ch, c some constant independent of h, O, O' . 

(5) Conservation: the reconstruction v should satisfy 

S{0)u'^ = [ v{x)dVix). 
Jo 
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For every remapping step at some time t, the old and new cells generate a layer of 
faces. For every cell TV that meets a cell O' with d{0' , O) < Ch, add a face O—^N 
and set 

Fo-.N^{ v{x)dV{x). 

J N 

By the assumptions, F has bounded stencil, is consistent, uniformly continuous and 
uniformly bounded; defining F^-tO —Fq^n renders it conservative. Moreover, 
it is clear that if the old and new grid satisfy the requirements outlined in the 
introduction and if the remapping steps are at least Q{h) apart, then the resulting 
space-time grid is quasiuniform. 

This technique will be discussed in more detail in future work (see jEllbj ) . 

5.4. Selfsimilar flow. Selfsimilar solutions, i.e. those that satisfy u{t, x) = u{st, sx) 
for all s > 0, arise in many important circumstances, such as shocks, contacts 
or rarefaction waves Riemann problems, or in |E1103j . A selfsimilar solution to 
ut + div/(it) = satisfies the system 

where £, = j are called similarity coordinates. To find asymptotically stable self- 
similar solutions, one can solve 

Ur + div^(/(M) — ^u) = —du 

(where r is a time- marching "pseudo-time" without physical significance). Adding 
the source term —du to numerical schemes in a conservative and consistent way (as 
defined in Section !?^ is easy; achieving stability is not difficult either. Theorem^] 
states that such a scheme delivers entropy solutions as long as it converges. 
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